nyc_border <- st_read("https://data.cityofnewyork.us/resource/wh2p-dxnf.geojson")
inundation_border <- st_read("https://data.cityofnewyork.us/resource/5xsi-dfpx.geojson")
# create bounding box to focus on the beach
riis_coords <- c(-73.893936,40.557853,-73.857586,40.578100)
names(riis_coords) = c("xmin", "ymin", "xmax", "ymax")
riis_bbox <- st_bbox(riis_coords)
riis_polygon <- st_as_sf(st_as_sfc(riis_bbox), crs = 4326)
# add topography limited to borders of beach and inundation zone
create_topo <- function(st, bbox, z_value) {
elev <- get_elev_raster(st, z = z_value, src = "gl1") # requires OpenTopography API key (https://opentopography.org/developers)
elev <- crop(elev, bbox)
elev_masked <- mask(elev, bbox) # remove NAs to reduce matrix dims
st_matrix <- raster_to_matrix(elev_masked)
return(st_matrix)
}
# higher zoom value = higher resolution
riis_matrix <- create_topo(nyc_border, riis_polygon, 14)
inundation_matrix <- create_topo(inundation_border, riis_polygon, 14)
zscale <- 10
n_frames <- 75
inundation_waterdepths <- max(inundation_matrix, na.rm = TRUE) - max(inundation_matrix, na.rm = TRUE) * cos(seq(0,2*pi,length.out = n_frames))
img_frames <- paste0("surge", seq_len(n_frames), ".png")
for (i in seq_len(n_frames)) {
message(paste(" - image", i, "of", n_frames))
riis_matrix |>
sphere_shade(texture = "imhof3", zscale = zscale) |>
add_water(detect_water(riis_matrix, cutoff = 1.3)) |>
add_shadow(ray_shade(riis_matrix, zscale = zscale)) |>
plot_3d(heightmap = riis_matrix,
zscale = zscale,
solid = FALSE,
shadow = FALSE,
water = TRUE,
watercolor = "#3289a0",
theta = -5,
phi = 65
)
render_water(inundation_matrix, waterdepth = inundation_waterdepths[i]/zscale,
watercolor ="#3289a0", zscale = zscale, remove_water = FALSE)
render_snapshot(img_frames[i])
rgl::clear3d()
}
magick::image_write_gif(magick::image_read(img_frames),
path = "riis_inundation.gif",
delay = 5/n_frames)
file.remove(list.files(pattern = "surge"))